Motivation of the study

This study aims to calibrate Leaf Area Density (LAD), which is defined as the area of leaves per unit of tree crown volume (in \(m^2.m^{-3}\)). This parameter forms the basis of numerous process-based forest models, such as physiological models that consider photosynthesis and transpiration (e.g. SurEau by Cochard et al., 2021), and forest gap models that estimate light interception (e.g. ForCEEPS by Morin et al., 2021, and Samsara2 by Courbaud et al., 2015).

Although the leaf area density can be directly estimated from a sample of trees (Stadt & Lieffers, 2000), this method is laborious and difficult to replicate for every tree in a stand. Therefore, indirect methods can be employed, such as model inversion using a hemispherical photograph (Courbaud et al., 2003) or the more recent LiDAR methodology (Wei et al., 2020). However, there is a significant lack of research on LAD in the scientific literature, as there are no large-scale databases containing values for a wide range of European species. Furthermore, the few reported species-specific values vary considerably between studies (Ligot et al., 2014).

The aim of this study is to estimate leaf area density (LAD) for nine different species/groups: Fagus sylvatica, Picea abies, Quercus sp., Carpinus betulus, Pseudotsuga menziesii, Pinus sylvestris, Abies alba, Larix decidua and ‘other species’. Additionally, Nock et al. (2008) demonstrated intraspecific variation in LAD, with intracrown leaf area decreasing by up to 40% with increasing tree age. This led us to consider the effect of diameter as well. Finally, as crown dimensions are very sensitive to local competition (Touzot et al., 2025), we hypothesise that LAD may also vary in response to competition.

Methodology

In this study, we used a model inversion methodology, using the ray-tracing model SamsaraLight (Courbaud et al. 2003). We want to estimate the SamsaraLight tree-level parameter LAD using measurements of transmitted light \(PACL_{obs}\) (for percentage of above light canopy), based on hemispherical photographs, within 45 plots of different composition and structures. For each plot, trees are inventoried with their crown dimensions, allowing to derive an estimated transmitted light \(PACL_{est}\) on virtual sensors with the ray-tracing model SamsaraLight. By doing so, for each virtual sensor of each plot, we have the observed transmitted light \(PACL_{obs}\), and the estimated transmitted light given a tree-level value of LAD \(PACL_{est}\).

The aim of the calibration is to find a tree-level LAD model that minimize the error between \(PACL_{obs}\) and \(PACL_{est}\). To do so, we defined the LAD for a tree \(i\), of species \(s\) in a plot \(p\) as \(LAD_{i,s,p} = a_{0,s} + a_{0,p} + a_{1,s}.DBH_i + a_{2,s}.BAH_i + \epsilon\), where \(a_{0,s}\) is a species-specific intercept, \(a_{0,p}\) is a random hierarchical site effect, nested within an origin effect (4 origins of the plots within UCLouvain, INRAe, Cloture and IRRES), \(DBH_i\) is the diameter at breast height of tree \(i\) , with the associated estimated species-specific coefficient \(a_{1,s}\), \(BAH_i\) the competition index of the tree \(i\) (i.e. the sum of basal area of trees higher than the tree \(i\)), with the associated estimated species-specific coefficient \(a_{2,s}\), and \(\epsilon\) the normal error term. We used a Bayesian approach to estimate the posterior distribution of parameters that minimize the log-likelihood between observed and simulated PACL. As computation of \(PACL_{est}\) needs an external call of SamsaraLight model, we used the R package ‘BayesianTools’ to run a MCMC sampling of parameters with the “DREAMzs” sampling algorithm.

For this study, we used the SamsaraLight model with the R package SamsaraLight (Beauchamp et al., work in progress, intended for publication in the Journal of Open Source Software). The SamsaraLight model was initially developed in Java and included in the Capsis simulation platform (SamsaraLightLoader by Gauthier Ligot). It could be called from R software via the “RCapsis” R package (Fortin et al., 2021). To speed up calculation times and make SamsaraLight easier to use in an R environment, we created the “SamsaraLight” R package. The source code is written in C++ to enable fast calculations and communication with the R environment is handled by the Rcpp package (Eddelbuettel & Balamuta, 2018). The “SamsaRaLight” R package reduces calculation time compared to the previous method (the RCapsis R package calling the Java model SamsaRaLightLoader) as it does not require the creation of a Java server nor Java objects. It also facilitates the use of SamsaraLight within R, as it does not require the creation of an inventory file, instead using R objects and functions directly. The package is functional and stored in a GitHub repository, but it is still subject to private restrictions as discussions regarding co-authorship and the intellectual property rights of the source code need to be concluded.

Calibration plots

The 45 calibration plots comes from 4 different sources we wiil call IRRES (9 plots, Gauthier Ligot), CLOTURE (23 plots, Gauthier Ligot), INRAe (10 plots, Philippe Balandier) and UCLouvain (3 plots, Mathieu Jonard).

IRRES

CLOTURE

INRAe

UCLouvain

Preliminary analysis

As a preliminary analysis, we tried to fit a mean LAD value for each sensor of each site. To do so, we run the SamsaraLight model on each stand by fixing the same crown LAD values of each tree, and for each sensor, we computed the residuals between the estimated PACL value on the virtual sensor, and the observed PACL on the field for this sensor. We did so for 500 values of LAD from 0.01 to 5. By doing so, we could find for each sensor the LAD that minimizes the residuals.

Many sensors did not converged (i.e. best fitting LAD was 5 as the residuals function showed an asymptotic form). Indeed, we could increase as much as we want the trees’s LAD, the estimated PACL could not decrease as much as the observed one. PACL is principally linked to the crown geometries, and secondary by the crown LAD. Thus, there are two cases where crown LAD did not influence PACL: (i) for open stands, where the high light on the ground comes principally from unobstructed rays, (ii) for highly crowded stands, where the very low light on the ground comes from the high number of intercepted crowns. Another reason could be linked to the wrong representation of leaves and crown above the sensor, where a simple branch can obstruct the sensor in the field where we considered a simple crown shape that result in incapabilities of representing such a light obstruction with the SamsaraLight model. Consequently, we remove the sensors that did not converged (465/1121 sensors), leading to remove consequently 3 sites where all sensors did not converged: Cloture11, Cloture15, Cloture2.

The mean LAD over all sensors that converged was about 0.94

Bayesian calibration

Observe the effect of sensors filtering

A first Bayesian approach would be to also fit a mean LAD value on the sensors that converged with the previous method, by considering a hierarchical random effect (site into origin, mandatory according to the structure our our dataset).

The model can be written as: \(LAD_{p} = a_{0} + a_{0,p} + \epsilon\), where \(a_{0}\) is the LAD intercept and \(a_{0,p}\) is a random hierarchical site effect, nested within an origin effect (4 origins of the plots within UCLouvain, INRAe, Cloture and IRRES).

Computation time

##   id_mod n_chains n_iterations filter_sensors time_hours
## 1      1        1        50000          FALSE   44.84417
## 2      2        1        50000           TRUE   35.34621

Compare models and intercept predictions

##   id_simu filter_sensors lad_qt025 lad_median lad_qt975
## 1     1_1          FALSE 1.2158727  1.4280842 1.5898866
## 2     1_2           TRUE 0.7916838  0.8749361 0.9477098
##   id_simu filter_sensors  mae_qt025 mae_median  mae_qt975 rmse_qt025
## 1     1_1          FALSE 0.06012604 0.06096777 0.06186939 0.08353675
## 2     1_2           TRUE 0.03579096 0.03637210 0.03720016 0.04776576
##   rmse_median rmse_qt975
## 1  0.08422148 0.08481025
## 2  0.04815073 0.04878985

Observe convergence

Is there an effect of species and tree diameter ?

##                      1_2                      2_1                      3_1 
##                "control"                    "dbh"                "species" 
##                      3_2                      4_1                      4_2 
##          "dbh + species"            "dbh*species" "dbh*species with covar"

Computation time

##   id_simu n_chains n_iterations intercept_per_sp dbh_effect dbh_per_sp
## 1     1_2        1        50000            FALSE      FALSE      FALSE
## 2     2_1        1        50000            FALSE       TRUE      FALSE
## 3     3_1        1        50000             TRUE      FALSE      FALSE
## 4     3_2        1        50000             TRUE       TRUE      FALSE
## 5     4_1        1        50000             TRUE       TRUE       TRUE
## 6     4_2        1        50000             TRUE       TRUE       TRUE
##   consider_covariance time_hours
## 1               FALSE   35.34621
## 2               FALSE   35.34444
## 3               FALSE   35.32820
## 4               FALSE   35.33053
## 5               FALSE   35.38885
## 6                TRUE   35.42728

Compare models

##   id_simu                str_mod intercept_per_sp dbh_effect dbh_per_sp
## 1     1_2                control            FALSE      FALSE      FALSE
## 2     2_1                    dbh            FALSE       TRUE      FALSE
## 3     3_1                species             TRUE      FALSE      FALSE
## 4     3_2          dbh + species             TRUE       TRUE      FALSE
## 5     4_1            dbh*species             TRUE       TRUE       TRUE
## 6     4_2 dbh*species with covar             TRUE       TRUE       TRUE
##   consider_covariance  mae_qt025 mae_median  mae_qt975 rmse_qt025 rmse_median
## 1               FALSE 0.03579096 0.03637210 0.03720016 0.04776576  0.04815073
## 2               FALSE 0.03651481 0.03714343 0.03828323 0.04807089  0.04896405
## 3               FALSE 0.03608892 0.03742451 0.03863126 0.04805010  0.04880111
## 4               FALSE 0.03641156 0.03681149 0.03913553 0.04787540  0.04825433
## 5               FALSE 0.03631365 0.03712481 0.03789914 0.04783215  0.04855410
## 6                TRUE 0.03616986 0.03693406 0.03769571 0.04774473  0.04823613
##   rmse_qt975               looic
## 1 0.04878985 -2085.44 (SE 44.64)
## 2 0.05028210 -2062.23 (SE 44.91)
## 3 0.04967238 -2056.22 (SE 43.22)
## 4 0.05031939  -2065.8 (SE 42.55)
## 5 0.04959985 -2063.38 (SE 42.07)
## 6 0.04888938  -2082.4 (SE 43.52)

Observe convergence

LAD predictions

For a mean species at a mean DBH

##                  str_mod lad_qt025 lad_median lad_qt975
## 1 dbh*species with covar 0.6510269  0.7324680 0.8259533
## 2                    dbh 0.6209563  0.7732074 0.8311620
## 3          dbh + species 0.6390058  0.7918842 0.8643817
## 4            dbh*species 0.7287725  0.8423970 0.8835583
## 5                species 0.7749527  0.8614020 1.1353462
## 6                control 0.9356714  0.9961420 1.1183788

Effect of DBH on mean LAD

Species-specific at a mean DBH

##                   str_mod               species lad_qt025 lad_median lad_qt975
## 1                 control            Abies_alba 0.9356714  0.9961420 1.1183788
## 2                 control      Carpinus_betulus 0.9356714  0.9961420 1.1183788
## 3                 control       Fagus_sylvatica 0.9356714  0.9961420 1.1183788
## 4                 control         Larix_decidua 0.9356714  0.9961420 1.1183788
## 5                 control           Picea_abies 0.9356714  0.9961420 1.1183788
## 6                 control      Pinus_sylvestris 0.9356714  0.9961420 1.1183788
## 7                 control Pseudotsuga_menziesii 0.9356714  0.9961420 1.1183788
## 8                 control            Quercus_sp 0.9356714  0.9961420 1.1183788
## 9                 control                 other 0.9356714  0.9961420 1.1183788
## 10                    dbh            Abies_alba 0.6209563  0.7732074 0.8311620
## 11                    dbh      Carpinus_betulus 0.6209563  0.7732074 0.8311620
## 12                    dbh       Fagus_sylvatica 0.6209563  0.7732074 0.8311620
## 13                    dbh         Larix_decidua 0.6209563  0.7732074 0.8311620
## 14                    dbh           Picea_abies 0.6209563  0.7732074 0.8311620
## 15                    dbh      Pinus_sylvestris 0.6209563  0.7732074 0.8311620
## 16                    dbh Pseudotsuga_menziesii 0.6209563  0.7732074 0.8311620
## 17                    dbh            Quercus_sp 0.6209563  0.7732074 0.8311620
## 18                    dbh                 other 0.6209563  0.7732074 0.8311620
## 19          dbh + species      Carpinus_betulus 0.5023566  0.6446167 0.6735185
## 20          dbh + species                 other 0.5614272  0.6637817 0.7456580
## 21          dbh + species Pseudotsuga_menziesii 0.6102095  0.7110653 0.7876103
## 22          dbh + species         Larix_decidua 0.6335558  0.7189264 0.8260378
## 23          dbh + species           Picea_abies 0.6218570  0.7635397 0.8355841
## 24          dbh + species       Fagus_sylvatica 0.6454982  0.7791721 0.8664678
## 25          dbh + species            Quercus_sp 0.6785826  0.7830967 0.8685027
## 26          dbh + species      Pinus_sylvestris 0.6898058  0.8031001 0.8785450
## 27          dbh + species            Abies_alba 0.6788627  0.8835034 1.0595781
## 28            dbh*species Pseudotsuga_menziesii 0.2485093  0.5269367 0.6388746
## 29            dbh*species                 other 0.4348611  0.6126360 0.6871356
## 30            dbh*species      Carpinus_betulus 0.2968384  0.6131672 0.8969416
## 31            dbh*species       Fagus_sylvatica 0.7204264  0.7588783 0.9182798
## 32            dbh*species            Abies_alba 0.7529801  0.8162976 0.9386565
## 33            dbh*species            Quercus_sp 0.6338826  0.8935212 1.1547369
## 34            dbh*species      Pinus_sylvestris 0.6995661  0.8947829 1.0275327
## 35            dbh*species         Larix_decidua 0.7392963  0.9683830 1.1205625
## 36            dbh*species           Picea_abies 0.8665411  1.2533810 1.4548950
## 37 dbh*species with covar      Carpinus_betulus 0.5210783  0.5809099 0.7419354
## 38 dbh*species with covar                 other 0.4464262  0.5994406 0.6806758
## 39 dbh*species with covar       Fagus_sylvatica 0.5707477  0.6682424 0.7768807
## 40 dbh*species with covar         Larix_decidua 0.6551333  0.7443774 0.8234916
## 41 dbh*species with covar            Abies_alba 0.7099401  0.7737799 1.1390419
## 42 dbh*species with covar Pseudotsuga_menziesii 0.6841458  0.7889470 1.2581770
## 43 dbh*species with covar            Quercus_sp 0.6755639  0.8325104 0.8917199
## 44 dbh*species with covar      Pinus_sylvestris 0.7395437  0.8910407 1.0222246
## 45 dbh*species with covar           Picea_abies 0.7740789  0.9349503 1.1605455
## 46                species                 other 0.3878365  0.6418581 0.7557415
## 47                species         Larix_decidua 0.5912653  0.6759419 0.7858099
## 48                species       Fagus_sylvatica 0.6413617  0.7154333 0.8369031
## 49                species Pseudotsuga_menziesii 0.6353227  0.7440858 0.8227895
## 50                species            Quercus_sp 0.7105984  0.7558197 0.9065006
## 51                species           Picea_abies 0.6394387  0.7795805 0.8902374
## 52                species            Abies_alba 0.7586773  0.8153506 1.1022361
## 53                species      Pinus_sylvestris 0.7943207  0.8653259 1.0237139
## 54                species      Carpinus_betulus 0.8616282  1.0184473 2.0344043

Effect of DBH for each species

## # A tibble: 1 × 3
##   rho_lower rho_median rho_upper
##       <dbl>      <dbl>     <dbl>
## 1    -0.706     -0.547    -0.331